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A new class of artificial atoms, such as synthetic nanocrystals or vortices 
in superconductors, naturally self-assemble into ordered arrays. This prop- 
erty makes them applicable to the design of novel solids, and devices whose 
properties often depend on the response of such assemblies to the action of 
external forces. Here we study the transport properties of a vortex array 
in the Corbino disk geometry by numerical simulations. In response to an 
injected current in the superconductor, the global resistance associated to 
vortex motion exhibits sharp jumps at two threshold current values. The 
first corresponds to a tearing transition from rigid rotation to plastic flow, 
due to the reiterative nucleation around the disk center of neutral dislocation 
pairs that unbind and glide across the entire disk. After the second jump, 
we observe a smoother plastic phase proceeding from the coherent glide of a 
larger number of dislocations arranged into radial grain boundaries. 



The production of ordered self-assembled structures of various materials as 
diverse as synthetic nanocrystals, magnetic colloids, charged particles in 
Coulomb crystals, proteins and surfactants, or vortices in type II supercon- 
ductors and in Bose-Einstein condensates, has attracted much interest for 
various fundamental and practical reasons which are ultimately concerned 
with their collective properties (optical, magnetic, mechanical, or transport 
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properties) PQI21E1- in particular, much experimental and theoretical effort 
has been devoted to characterizing the phase diagram of type II supercon- 
ductors 4j. Depending on the value of the magnetic field H, temperature 
T, and sample preparation, vortices can either form a crystal jS], which at 
higher temperatures melts into a liquid El M > or ; due t° quenched disor- 
der, they can be found in more complex phases, such as the vortex glass [TU] , 
the Bose glass [TT] or the Bragg glass [121 H3] ■ Of special importance is the 
non-equilibrium response of vortex matter to the flow of an external current 
[EJUS], because the dissipative motion of the vortices induces an undesirable 
macroscopic resistance. The moving phase can be as simple as the collective 
motion of an elastically deforming vortex crystal, or can be more complex, 
such as in plastic [10] or in channel vortex flow [T7| . 

Transport experiments in superconductors are often performed in a strip 
geometry: current is injected in one side and removed from the opposite 
side [inillZj. The current- voltage (I-V) curve provides an indirect measure of 
vortex dynamics, because vortex motion induces an electric field proportional 
to the vortex velocities. It has been noticed that the sample boundary has an 
important effect on the moving vortex phase complicating the interpretation 
of the results jTHj II 9j : a problem that is overridden in the Corbino disk 
geometry ^3 E] • There the current is applied at the disk center 

and flows radially towards the boundary. Vortices tend to move in concentric 
circles without crossing the sample boundaries, avoiding edge contamination. 

Vortices in the Corbino geometry experience a force gradient ^5], and 
thus exhibit intriguing dynamic phases as a function of T, H, and /. The 
vortex velocity profiles have been evaluated after measuring the voltage drop 
across a series of contacts placed radially on a YBa2Cus07_5 disk |21j . For 
low currents and temperatures, all the vortices were found to move as a 
rigid solid, giving rise to a linear velocity profile v(r) = Qr, where r is the 
distance from the disk center. Above a threshold current Jo, the vortex crys- 
tal cannot sustain the shear stress induced by the resulting inhomogeneous 
Lorentz force, and the response becomes plastic. Finally, above the vortex 
lattice melting temperature T M , the velocity profile is fluid-like and decays 
as v(r) ~ 1/r. A theoretical model of vortex flow in the Corbino geometry 
has been analyzed, and shows the shear yielding as a dislocation unbinding 
transition [221 123] ■ Plastic flow would appear as soon as the current-induced 
shear stress is large enough to separate an existing pair of bound dislocations. 

In this paper, we first study transport in the Corbino disk by T = 
molecular dynamics (MD) simulations of interacting vortices [2H I2SH2E1- As 
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in the experiments, for low currents we find a linear velocity profile that cor- 
responds to the rigid rotation of the vortex lattice. Above a threshold current 
Iq, the profile ceases to be linear, indicating the onset of plastic flow. Our 
simulations enable a close inspection of the lattice topology, which unveils 
the microscopic origin of this dynamic transition: at and above Iq, new dis- 
location pairs are created, mainly within the highly strained central region, 
which readily unbind and glide along all possible crystalline directions giving 
rise to plastic flow. These processes occur repeatedly, yielding a strongly 
fluctuating voltage noise, which is reminiscent of the intermittent behavior 
observed in plastically deforming crystals j2l2Zj- For currents larger than a 
second threshold Ii, we observe that voltage fluctuations decrease and vor- 
tices end up moving in uncorrelated annular channels, displaying a laminar 
1/r velocity profile. In this regime, we find a larger amount of dislocations 
in the crystal, most of them forming radial grain boundaries that span the 
entire disk and glide in the tangential direction. The crystal reorientations 
associated to the presence of these grain boundaries make possible a steady 
regime of plastic deformation in the azimuthal direction. In addition, the 
exponential screening of shear stress produced by grain boundaries |2Bj en- 
hances the number of nucleation events. The number of dislocations after 
the second jump at I\ reaches a maximum value, corresponding to the pres- 
ence of quite densely packed grain boundaries. We observe that the plastic 
deformation of the crystalline film proceeds with the glide motion of these 
grain boundaries. 

We measure the plastic threshold current Iq for different values of the 
disk radius D and vortex number N, and later show (in Fig. EJ) that the data 
collapse into a single curve when plotted against the crystal average lattice 
spacing a. As predicted previously [231^0 is proportional to the vortex crystal 
shear modulus c 66 . Indeed, the curve follows very closely the dependence of 
C66 on the lattice spacing or field strength A simple evaluation of the 
energy cost of a new dislocation pair provides a good quantitative estimate 
of this threshold current. 

In the Corbino geometry, a disk-shaped superconductor is placed in a 
magnetic field parallel to the disk axis, and a current I is injected at a metal 
contact in the disk center and removed at the disk boundary. Thus, the 
current density inside the disk is given by J(r) = rl/(2nrh), where h is the 
thickness of the specimen. This radial current generates an azimuthal Lorentz 
force acting on the vortices fifr) — @&oJ{ r )/c, where $o is the quantized 
flux carried by the vortices, c is the speed of light, and 9 is the azimuthal 
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versor. In addition, a pair of vortices interact with each other via a long- 
range force f OT (r) = AK"i(|r|/ A)f, where A = <E>q/(87t 2 A 3 ), A is the London 
penetration length, and K\ is a Bessel function |3U]. Taking into account 
these interactions, we solve numerically the dynamics of iV vortices confined 
in a disk of radius D (see the Methods section, describing the technical details 
of the simulations). 

As in previous experiments [21], we first measure azimuthal velocities v 
and compute the variation of the velocity profile v(r) with I (see Fig. Q). 
At low currents, the velocity profile follows a linear law, v(r) = lot with 
angular velocity Q oc I. This corresponds to a rigid rotation of the vortex 
crystal, as observed experimentally |2T] . Above a threshold Jo, the velocity 
profile starts to deform, with large tangential velocities in the center, which 
then decay towards the boundary. At higher currents, the profile becomes 
smoother, decaying as 1/r, characteristic of a laminar response. 

To better identify the transitions in the system rheology, we measure the 
variations of the flow resistance R = J2i v i/I with I. After an initial tran- 
sient, the resistance reaches a steady state, which fluctuates strongly in the 
plastic regime and is much smoother in the solid and laminar phases. Figure 
I2ja) shows the steady-state resistance for different values of N (at D — 18) 
as a function of I. We have normalized these curves by the corresponding 
number of vortices iV to better visualize their characteristic features. The 
curves show a first sharp jump around Iq corresponding to the breakdown 
of the linear velocity profile, and a smaller jump at ii, indicating the on- 
set of the hyperbolic profile. The final plateau scales with the number of 
moving vortices N. Indeed in this laminar regime, a scaling factor of N/D 
follows from a simple continuum approximation with a constant density of 
vortices. A direct inspection of the topology of the lattice allows the nature 
of the transitions to be clarified. We construct the Delaunay triangulation of 
the vortex positions in the disk to characterize their topology. Most of the 
vortices are sixfold, as in a perfect triangular lattice such as the Abrikosov 
lattice. A pair of fivefold and sevenfold neighboring vortices identifies an edge 
dislocation in the lattice, a topological defect characterized by its Burgers 
vector b [2E]- Dislocations produce long-range stress and strain fields in the 
host crystal, experience the so-called Peach-Koehler force due to the local 
stress, and move mainly by gliding along the direction of b |2E] ■ 

For I < Iq, all vortices within the bulk of the disk are six-fold whereas 
a large number of fivefold- and sevenfold-coordinated vortices are only ob- 
served along the boundary. These are geometrically necessary dislocations 
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and disclinations, which need to be present in order to adjust a triangular lat- 
tice into a circular geometry. Like the resistance, the number of five/sevenfold 
coordinated vortices fluctuates around a steady average value after an initial 
transient. Figure EJb) shows the behavior of the average steady number of 
fivefold vortices n 5 as a function of the current. We have subtracted the 
average number of geometrically necessary boundary defects n 5 (0), and nor- 
malized the curves by its maximum value n™ ax — n 5 (0) to better visualize 
their main features. The curves in Fig. EJb) closely resemble the behavior of 
the resistance in Fig. |2fa), with jumps at Iq and I\ (see also Supplementary 
Information, Fig. SI corresponding to a disk of radius D = 36). As the 
current overcomes Jo, new defects start to nucleate near the center of the 
lattice. Typically, we observe the reiterative formation of new dislocation 
dipoles (two dislocations with opposite Burgers vectors), that unbind and 
glide along the direction of their Burgers vector, in most cases towards the 
disk boundary (see Fig. Efa) and supplementary movies 1 and 2). To accom- 
modate the shear stress generated by the external current, the crystal should 
nucleate dislocations that are able to glide either radially or tangentially. 
Nevertheless, in the undistorted triangular lattice (or when the concentra- 
tion of free dislocations is low), the dislocations that are nucleated are the 
most elementary, with Burgers vectors along the three basic crystalline direc- 
tions. Of those, only the ones with b almost parallel to the radial direction 
can easily glide over long distances due to the Peach-Koehler forces involved. 

As time goes on, the dislocation flow process exhibits an erratic character, 
because dislocation pairs at short distances may annihilate each other or react 
to form a new dislocation; they assist the nucleation process at intermediate 
distances, and may even form various metastable structures. This intricate 
process is reflected by the strongly fluctuating flow resistance and the non- 
linear velocity profile. Figure 0] shows that at the onset of the plastic phase, 
the resistance noise power spectrum S(u) (where u is the frequency) displays 
a non-trivial power law decay lo~^ with (3 = 1.5. This is also apparent from 
the inset of Fig. which displays the relative standard deviation of the 
resistance fluctuations as a function of current. We can then conclude that 
in the Corbino disk, strong voltage noise is a fingerprint of the onset of plastic 
deformation in the lattice. 

Further increase of the current beyond I — I\ results in a smoother 
flow, with most vortices moving in uncorrelated concentric trajectories. Here 
dislocation flow appears to be quite peculiar because a significant number 
are settled in walls (grain boundaries) oriented along the radial direction, 
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which tend to glide coherently in the azimuthal direction. Grain boundaries 
produce the necessary deformations of the crystal that allow a stationary 
tangential flow (see Fig. 0{b) and supplementary movies 3 and 4). Besides, 
the long-range stress field generated by free dislocations is screened out with 
their arrangement into grain boundaries. The screening length grows with 
the average distance separating contiguous dislocations in the wall, that is, 
the grain boundary spacing. In a close-packed grain boundary, the spacing 
is of the order of a few crystal lattice constants a, and consequently, the 
number of dislocations in a radial grain boundary of size D grows as D/a. 
We indeed observe that the asymptotic number of dislocations in the laminar 
regime grows as n™ ax ~ D/a ~ y/N (see Supplementary Information, Fig. 
S2). We also observe that relative fluctuations of dislocation number and 
resistance (see inset in Fig. 0} are consequently reduced in this regime. Thus 
at and above I\, the rate of nucleation of new dislocations is big enough 
to ensure the formation and maintenance of these type of grain boundaries 
whose cooperative motion marks the new regime of plastic deformation. 

Motivated by the observation that the onset of plasticity corresponds to 
the reiterative nucleation of topological defects around the center of the disk, 
we estimate Jo by computing the energy cost of a dislocation dipole [23] • Very 
much like in any other ordinary crystal, the energy cost of an edge dislocation 
in the vortex lattice is made up of two contributions: the core energy F c and 
the elastic energy F e . The elastic energy cost F e is proportional to the shear 
modulus of the crystal, and grows with the logarithm of the system size 
D [2E]- The order of magnitude of the core energy per unit length can be 
estimated using the same variational argument that has been proposed [3T] 
for a vortex line in the xy-model. It turns out that F c is proportional to 
c 66^ 2 /(4vr) for each dislocation, where Cqq is the local shear modulus of the 
vortex lattice. 

The elastic energy cost of a new dislocation dipole, such as the ones 
observed in the plastic phase, is independent of the system size and grows 
with the relative distance between the dislocations — which, right at the 
moment of the creation event, is of the order of the lattice spacing — as 
6 2 c 66 /(27r) \n(r d /r c ), where r c is a variational short distance cutoff. Because 
both rd and r c are of the order of the lattice spacing, the core energy of 
the new dislocation dipole 2F C is the leading energetic contribution for the 
proliferation of new pairs. The inhomogeneous elastic shear stress and strain 
induced by the external current in the Corbino disk geometry have been cal- 
culated analytically in Ref. |23j. According to those results, the elastic energy 
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stored in a region around the disk center of size R d — of the order of the spa- 
tial extent of a dislocation dipole — is roughly equal to (R d IB/47ich) 2 7r/ (4c 66 ), 
where B = <fton is the average magnetic induction and n the areal density 
of flux lines. This energy is released by the formation of new dislocation 
pairs. On balancing the core energy with the elastic energy provided by the 
external current one can estimate the transition current 

_ Ay/2chbc m _ 4V2chb(j) c 66 
°" R d B ~ i^(87rA) 2 c 66 ' U 

where c^e = B<po/ '(8tt\) 2 is the long wavelength shear modulus in the con- 
tinuum limit. In the units of current used in the simulations this transition 
current is equal to J = V2b/ (4irR d )c e6 /cQ 6 . 

The vortex lattice is often considered as a conventional continuum elastic 
medium characterized by its compressional en, shear cqq, and tilt C44 moduli, 
disregarding its discrete nature. In the limit of very long wavelength distor- 
tions, however, the elastic properties are governed by local elastic moduli, 
which strongly depend on the strength of the magnetic field. Analytical re- 
sults for the local compressional, shear, and tilt moduli of a discrete vortex 
lattice as a function of the lattice spacing a are provided elsewhere [29J. 

We compute systematically the variation of Iq with iV and D (see the 
inset of Fig. EJ). When plotted as a function of the lattice spacing a = 

\J2itD 2 / V3N, the scattered data can be collapsed into a single curve (Fig. 
main plot) which follows the theoretical curve c 66 /c 66 . When multiplied 
by an overall constant C = 0.039 ± 0.003 the curve provides a good fit of the 
data. This constant is also in reasonable good agreement with its theoretical 
estimate C = V2b/ (<inR d ) (where b = a and the dipole extent R d ~ 2 — 3a). 

In conclusion, we have shown that the onset of plastic flow in a supercon- 
ducting disk occurs in close correspondence with the nucleation and motion 
of new dislocations in the vortex lattice. Once formed, dislocations glide 
parallel to their Burgers vector, releasing the shear stress concentration. On 
increasing the nucleation rate, dislocations arrange themselves into densely 
packed radial grain boundaries that cross the entire disk and tend to glide in 
a cooperative manner. The plastic response is characterized by a non-linear 
I-V curve and by strong fluctuations. These results have been obtained for 
a two dimensional geometry, and they could change when the thickness of 
the film is large. Our approach should be relevant, however, for thin films of 
other self-assembled nanostructures subject to the action of shearing forces 
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[111210]; and possibly for granular media [32] • Obviously, the precise value 
of the driving force thresholds, signaling the onset of a nonlinear response 
and the plastic deformation of the sample, would depend on the appropri- 
ate physical parameters characterizing the interaction among the constituent 
elements in each case. 

Methods 

We consider a set of N rigid vortices confined in a disk of radius D. The 
equation of motion for each vortex i at position 

Ydvi/dt = £ f OT ( ri - r s ) + f L (n), (2) 

3 

where T is an effective viscosity. We choose as units of space and time 
A and to — TA/A respectively, and we measure the current I in units of 
$o/(2vrc/iAy4). The N vortices are confined inside the disk by the external 
magnetic field and the sample edge barrier that we model by imposing an 
extra normal force on the vortices of the form fg = — gexp [—(D — r)/ro]/rof, 
with ro = 0.1 A and g/A — 1. A similar force is also imposed at the inner 
wall close to the disk center (at r = r ), thus avoiding the singularity of the 
Lorentz force at r = 0. 

The coupled Eqs. (J2J) for i — 1, ...,N are integrated numerically with an 
adaptive step size fifth-order Runge-Kutta method with precision 10~ 6 . We 
do not truncate the range of the vortex-vortex interaction since this leads 
to spurious fluctuations caused by the force discontinuities. We study the 
response of the system as a function of the applied current for different values 
of N, ranging from iV = 332 to iV = 2064, and D (D = 18A, 36A, 72A). 
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Figure 1: The velocity profiles as a function of the applied current for N = 
1032 vortices in a disk of radius D = 18. For I < I — 0.036 the profile is 
linear, corresponding to a rigid rotation of the lattice. For I > I the profile 
deforms, indicating plastic flow. At high drives, the profile simply decays as 
1/r, where r is the distance from the disk center. 



12 



0.11 




Figure 2: The resistance and the number of fivefold vortices as a function of 
the current, (a) The normalized resistance R/N as a function of the applied 
current for N vortices in a disk of radius D = 18. (b) The excursions in the 
average steady number of fivefold coordinated vortices n$ for different N by 
plotting (n 5 - n 5 (0))/« a:r - n 5 (0)). 
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Figure 3: A series of snapshots illustrating dislocation dynamics. The system 
consists of iV = 1032 vortices in a dis¥of radius D = 18 under two different 
applied currents: (a) / = 0.04, and (b) / = 0.3. Fivefold coordinated vor- 
tices are highlighted with colored circles. Different colors correspond to four 
different time steps separated by an interval At = 100. Initially they are red, 
then violet, blue, and cyan. As a reference, on each plot we also show the 
Delaunay triangulation for the initial time step considered. Notice the radial 
motion of some defects in the nlastic nhase (a), and the tangential elide of 
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Figure 4: The resistance noise power spectrum for iV = 1032 vortices in a 
disk of radius D = 18 under an applied current / = 0.0375. There is a region 
of power law decay well fit by c<j~ L5 . The inset shows the relative resistance 
noise standard deviation as a function of the current for D = 18 and D = 36. 
The peaks correspond to the transitions form solid to plastic and from plastic 
to laminar. 
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Figure 5: The threshold current I . In the inset we report the value of I 
as a function of the number of vortices iV for different disk radii D. When 
plotted as a function of the lattice spacing a, the data collapse into a single 
curve, which follows the decay of c 66 (main plot). 
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